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1.0 INTRODUCTION 


1 . 1 SUMMARY 

This report describes the work performed on spectral 
analysis of data from the C-Band radars tracking GEOS-II 
and on the development of a data compaction method for 
the GEOS-II C-Band radar data. 

The purposes of the spectral analysis study were to 
determine the optimum data recording and sampling rates 
for C-Band radar data and to determine the optimum method 
of filtering and smoothing the data. The optimum data re- 
cording and sampling rate is defined as the rate which in- 
cludes an optimum compromise between serial correlation 
and the effects of frequency folding. The goal in develop- 
ment of a data compaction method was to reduce to a mini- 
mum the amount of data stored, while maintaining all of 
the statistical information content of the non- compacted 
data. 

A digital computer program for computing estimates of 
the power spectral density function of sampled data was 
used to perform the spectral analysis study. See Section 

1.3.1 and Section A-2 of Appendix A for details. The spec- 
tral analysis was performed on raw data, before any filter- 
ing or data correction operations were performed, so that 
the data analyzed contained signal and noise components un- 
modified by data processing. 

The first study performed, described in Section 2.1.1, 
determined the signal cut-off frequency. It was determined 
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that for range, azimuth, and elevation data on tracks of 
GEOS-II, the signal component of the data is concentrated 
in frequencies below .05 Hertz. The second study on spec- 
tral content of the data, described in Sections 2.1.2 and 
2.1.3, was conducted to define the spectral characteristics 
of the noise content of the data, at frequencies above .05 
Hertz. The results are summarized in Section 3.1.1. 

Based on the results of the spectral analysis studies, 
it was determined that low-pass linear digital filtering of 
the raw data would reduce the noise content of the data sig- 
nificantly while not degrading the signal content. A digi- 
tal computer program to implement this approach was used 
to conduct studies on the properties of filtered data. The 
approach and the computer program are discussed in Section 
1.3.2 and in Section A- 3 of Appendix A. 

The data studies conducted using the filtering and 
compaction program are discussed in Section 2.2. These 
studies consisted of orbital reductions, using both filtered 
and unfiltered data, and a comparison of the orbits and of 
the data residuals from the two sets of data. The major 
conclusions to be drawn from these studies are that filter- 
ing improves the precision of the data, as demonstrated by 
the RMS of the short arc reductions or by examining the resid- 
ual plots, and that filtering does not significantly affect 
the accuracy of the data, as demonstrated by comparison of 
the orbital elements and by data extrapolation. Details of 
these conclusions are given in Section 3.1.2. 
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Studies to determine the possible serial correlation 
of filtered and unfiltered data were conducted using a 
program to compute the auto-covariance function, as describ- 
ed in Section 2.3. Two sets of data with known abnormalit- 
ies were investigated, the ship test data study as described 
in Section 2.3.2 and the PMRMR1 data as described in Section 
2.3.3. These studies indicate that calculations of the 
auto-covariance function can be used to investigate abnormal 
data. The investigation of the auto-covariance function of 
residuals from short-arc orbits, using both filtered and 
unfiltered data, is reported in Section 2.3.4. The results 
indicate that serial correlation of short arc residuals is 
not dramatically affected by data filtering. 

In order to implement the digital studies described, 
computer programs designed for the purpose were developed. 
These programs are described, in general, in Appendix A, 
and in Reference 2, and the detailed theory used is given 
in Reference 1. 

1 . 2 STUDY GOALS 

The goals of this segment of the effort are: 

a) Evaluation and analysis of the spectral content 
of the radar data. 

b) Development, implementation and documentation of 

a data compaction and reformatting scheme suitable 
for supplying C-Band data to the National Space 
Science Data Center. 

1.2.1 Spectral Analysis 

The first goal of this effort was to determine the 
spectral content of the Wallops AN/FPQ-6 and AN/FPS-16 
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radar measurement channels. Determination of spectral 
content, or power spectral density estimation, is concern- 
ed with decomposition of the time domain signal from a 
particular channel into periodic components, as a function 
of frequency, and determining the average power in each 
component. The term "power" used in this context refers 
to the mean squared amplitude of the component, in analogy 
to the mean squared amplitude of a sinusoidal voltage 
passed through a unit resistance. Spectral analysis is 
based on the mathematical principle that a set of data 
points representing a function of time can be transformed 
into a set of data representing a function of frequency. 

The frequency domain representation of the data is called 
the spectrum , or the spectral representation, of the time 
domain representation. The process of transforming data 
from the time domain to the frequency domain, and then 
evaluating the spectral characteristics, is called spectral 
analysis . The process of transformation from one domain 
to another does not add or subtract any information from 
the data. The information content of the frequency domain 
data is the same as that of the time domain data. 

Frequency domain analysis is a valuable tool when we 
wish to consider the effects of an additive measurement 
error. The radar data we are considering is the result of 
the measurement of a continuous physical process, and the 
measurements have had random errors added to them in the 
measurement process. In frequency domain terminology, 
which was derived from communications engineering, the 
error-free measurements are referred to as the signal and 
the additive errors are called noise . Each time domain 
point contains both the signal and noise components. In 
the frequency domain, however, the signal component tends 
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to be concentrated in the low frequency portion of the 
data, while the noise is spread over a wider range of the 
frequency band. To the extent that this is true, it is 
possible to distinguish between signal and noise in the 
frequency domain and to define methods for eliminating 
the noise. 

Spectral analysis is worth doing when the frequency 
domain representation can tell us more about the data than 
we can obtain from the time domain function. This con- 
dition is true when we are dealing with data which has been 
passed through a time invariant linear system. The effect 
of a time invariant linear system on each frequency component 
of the spectral representation of a function is multipli- 
cation by a constant phase angle. Since a radar instru- 
ment can be approximated by a time invariant linear system, 
spectral analysis of radar data can help in analyzing the 
effects of the radar on the data, and help to predict the 
effects of further data processing. 

Another piece of information to be obtained from, spec- 
tral analysis involves the tradeoffs in the choice of a 
sampling interval for data. The power spectrum of a data 
array and autocorrelation function of the data are a Fourier 
transform pair. That is, the power spectrum is the fre- 
quency domain , representation of the autocorrelation function. 
If the power spectrum of the noise has a cutoff at a fre- 
quency f c , the sampling theorem says that data sampled at 
an interval of more than l/(2f c ) is subject to aliasing or 
frequency folding . For a very narrow band signal, this 
aliasing does not degrade the data if the sampling interval 
is less than l/f c * This is a frequency domain constraint 
on the sampling interval, The time domain constraint 
on the sampling interval is based on the autocorrelation 
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function. Many methods of statistical data reduction, 
including orbit determination methods, contain error esti- 
mates based on an assumption of zero (or negligible) ser- 
ial correlation between data points. On this basis, the 
sampling interval should be the time for which the auto- 
correlation is sufficiently small. 

1.2.2 Data Compaction 

For much geodetic work, such as the estimation of 
geopotential coefficients or the recovery of station po- 
sitions, long arcs of data must be processed. At the nor- 
mal sampling rate for most radars of ten observations per 
second an excessive number of observations must be processed. 
The solution is not simply to discard data. A method 
based on simply selecting each Nth point and discarding 
N-l of every N points will throw away information with 
the discarded points. In the frequency domain represen- 
tation, this method of data compaction introduced aliasing 
or frequency foldover of the higher frequency noise com*- 
ponents of the spectrum. For example, if data containing 
a 1 hz noise component is sampled at 1/second, the 1 hz 
component becomes a bias error. 

A data smoothing or filtering scheme must be used 
so that the frequency of observations, or sampling rate, 
is significantly reduced while the majority of the orig- 
inal information is retained. The design goal is to re- 
duce to a minumum the amount of data that is necessary 
to maintain the statistical information content of the 
non-compacted data. The filtering operation must decrease 
the bandwidth of the noise content of the data, so that 
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the cut-off frequency f is reduced. Then the data sel- 
ection interval can be chosen to avoid frequency foldover 
effects and to minimize serial correlation of the data. 


1.3 TECHNICAL APPROACH 

1.3.1 Power Spectrum Estimation by Fourier Methods 

The power spectrum is an idealized mathematical 
concept and in the strictest sense is defined only for 
an infinite set of data. From any given finite set of 
data we can only obtain an estimate of the power spec- 
trum of a hypothetical set of infinite data from which 
the finite set is a sample. 

An intuitive, non-mathematical , definition of the 
power spectrum is based on the concept of passing the 
function through a parallel array of very narrow band 
pass filters. Each filter is an idealized perfect fil- 
ter, with unit gain and zero phase shift in a narrow 
band, and zero gain outside that band. The pass bands 
of the set of filters are non-overlapping, and cover the 
frequency band from zero to infinity. Then the output 
of each filter is squared, and the average power over 
time is obtained by integrating and obtaining a time av- 
erage. If the band width of each filter is allowed to 
approach zero, and the time over which the power is 
averaged approaches infinity, the resulting output is 
the power spectrum of the input function. 

Two approaches have been used in the past to ob- 
tain spectral density estimates from sampled data. The 
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direct approach uses the digital Fourier transform to 
obtain a spectral representation of the data, then 
obtains the estimate of power spectral density by squaring 
the spectrum and averaging to smooth the estimated power 
spectrum. Various methods of smoothing of the power 
spectrum are used. Smoothing or averaging of the spectrum 
is necessary because the unsmoothed squared spectral 
function, or periodogram , is a statistically unstable es- 
timator of the spectral density. 

The indirect method of obtaining an estimated power 
spectrum is the method most widely described in the 
literature, because it was the most widely used until a 
few years ago. The indirect method consists of obtain- 
ing an estimate of the auto- covariance function of the 
data by taking mean lagged products, then Fourier trans- 
forming the estimated auto-covariance function to obtain 
an estimated power spectral density function. The in- 
direct method is based on a theorem of Fourier analysis 
which states that the auto-covariance function and the 
power spectral density are a Fourier transform pair. 

The reason for the widespread use of the indirect 
method is not that it is theoretically better, because 
the indirect and the direct methods are theoretically 
equivalent. The indirect method was used because, before 
discovery of the fast Fourier transform, the direct method 
involved much too much computing time. Since the use of the 
fast Fourier transform has become widespread, the indirect 
method is being used less. 

For this project, the direct method was used. The 
calculation of the periodogram is implemented by use of 
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the fast Fourier transorm method. Two types of smoothing 
of the periodogram were programmed and used, a linear smoo- 
thing to increase resolution on the signal components of 
the data and a quadratic smoothing to decrease the var- 
iance of the estimation of the noise components. The 
theory of the methods used is discussed in Reference 1. 


1.3.2 Data Compaction Using Digital Filtering 

The method selected for use in preparing input data 
for data compaction was the use of linear digital filters. 
In the following discussion, the term "digital filter" 
rather than "linear digital filter" will often be used, 
but in general, the reference is to linear digital fil- 
ters . 


This distinction is made because any mathematical 
operation which changes the frequency content of a set 
of data can be considered to be a digital filter. For 
example, if an orbit is computed using measured range, az- 
imuth and elevation data and then the (R, A, E) var- 
iables are computed from that orbit at the measurement 
times, those computed variables can be considered to 
be the result of filtering the original data by use of 
the orbital equations. 

The decision to use linear digital filters for the 
data compaction method, rather than use a more elaborate 
"filter" such as a short arc orbit or a state variable 
dynamic model, was based on a desire not to introduce 
a priori assumptions into the analysis at this early 
stage. Any sophisticated data filtering method such 
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as fitting to an orbit involves deciding on a system model 
which can in any case only approximate the actual dynamics 
of the system being observed. The approximations intro- 
duced into the model can produce trends and biases into 
the resulting data. For example, the effects of station 
survey errors and gravity model errors can introduce trends 
and biases into an orbit reduced from the data. An inves- 
tigator desiring to use the data at a later date would not 
be able to remove these trends and biases if he desired 
to use a different dynamic model. 

The use of a linear digital filter in data prepar- 
ation involves only a very simple model of the data which 
is valid over a much shorter time interval than an orbit- 
al pass. That model is contained in the assertion that the 
signal portion of the data is concentrated in the low fre- 
quency part of the spectrum, and that the high frequency 
part of the spectrum is all noise. Then a low pass fil- 
ter, which acts as an identity operator on the low fre- 
quency components and effectively removes all high fre- 
quency components, will remove much of the noise without 
changing the signal. 

A linear digital filter operates on an input sequence 
x to produce an output sequence y. If the previously 
values of y are used in the computation of the next value, 
the filter is called recursive . If only input sequence 
values are used in computing each value of the output se- 
quence, the filter is non-recurs ive . 

For the applications in which we will be interested, 
non-recursive filters have the most general use. Recursive 
filters are more generally applicable to real time appli- 
cations, in simulation or in digital control systems. 
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We are concerned here with a da-ta processing application, 
in which the data to be filtered is stored in an array, 
with no distinction between "past", "present" or "future". 
In this situation, non- recurs ive filters can be used and 
are more versatile and general in application than re- 
cursive filters. 

The specific use for the digital filters will be 
as low-pass filters, to eliminate the noise component 
of the data while retaining the signal component. The 
filters to be used must be highly selective, with a very 
sharp cut-off above a specified frequency. In order to 
attain this high selectivity with a very narrow band 
low-pass filter, a large number of data input points 
must be filtered to obtain each output point. The 
method used to implement these filters employs the Fourier 
Convolution theorem. The process of applying a non-recur- 
sive digital filter to a set of data is mathematically 
equivalent to obtaining the convolution of the filter coe- 
fficients with the set of data. The Fourier Convolution 
theorem states that the convolution of two functions is 
the inverse Fourier transform of the products of the 
Fourier transforms of the two functions. This method is 
used for digital filtering because the employment of the 
fast Fourier transform algorithm makes it much more effic- 
ient than the direct convolution. 

A disadvantage in the use of non-recursive filters 
is loss of data points at both ends of the data span. The 
non-recursive filter can also be regarded as a weighting 
function, consisting of a symmetric array of coefficients 
which are multiplied by the input data points to produce 
one output data point at the center of the array. It is 
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impossible to obtain output data at the ends of the input 
data array because the weighting function must be applied 
to the input data symmetrically about the output data point. 
The number of data points lost in the filtering process 
is one half of the number of points in the filter minus 
one. The disadvantage of loss of data points is increased 
because the number of coefficients in the filter is de- 
pendent on the shaprness of the required cut-off. To 
a close approximation, the minimum number of points in 
the filter is inversely proportional to the cut-off fre- 
quency desired. Thus a tradeoff between cut-off frequency 
and data points lost must be considered. Such a trade-off 
is described in Section 3.2. 


1.3.3 Use of Finite Fourier Transform 

The theory of the methods used for power spectrum 
estimation and for computation of digital filters is 
based on the digital Fourier transofrm, which has been 
known and used for many years. The actual numerical 
methods used are based on a special way of calculating 
the digital Fourier transform, which was developed only 
a few years ago and which is called the fast Fourier 
transform. The fast Fourier transform algorithm is only 
a very efficient method of calculating the digital Fourier 
transform, and does not add anything the the theory of what 
can be done. In practice, because of the great decrease 
of computing time brought about by its use, the fast 
Fourier transform actually makes possible the application 
of Fourier analysis methods to routine problems. 
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2.0 RADAR DATA STUDIES PERFORMED 


The computer programs described in Appendix A were used 
to perform studies on the C-Band radar data. Most of the 
radar data used as inputs to the programs was taken from 
the GEOS- 1 1 tracking records. However, since the GEOS- 1 1 
tracking data is almost all beacon track, data from other 
satellites was used for some of the studies of skin track 
characteristics . 


2.1 SPECTRAL ANALYSIS STUDIES 

The spectral analysis studies were carried out using 
the SPECDEN program, described in Appendix A. The goal of 
these studies was to determine the spectral characteristics 
of the radar measurement data from the Wallops AN/FPQ-6 
and AN/FPS-16 radars. All of the spectral analysis studies 
were carried out on data containing both signal and noise. 
The data was not filtered or modified in order to remove 
the signal component before analyzing the noise character- 
istics. The basis of the decision not to prefilter the 
data in any fashion, before spectral analysis, was the 
desire to avoid any possibility that the spectral content 
of the data had been modified by the filtering process. The 
investigation was divided into two phases, first, to define 
as well as possible the cut-off frequency of the signal 
component of the data, and second to obtain the spectral 
characteristics of the noise component of the data. 


2.1.1 Signal Cut-Off Determination 

The spectral analysis of any data containing both 
signal and noise must be planned on the basis of some know- 
ledge of the characteristics of the data. This is a 
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seeming paradox: that the analysis cannot proceed until 

the results of the analysis are known. The paradox must 
be resolved by performing some preliminary analysis based 
on general knowledge of the data to be analyzed. The 
signal cut-off determination studies reported here were', 
part of that preliminary analysis. 

It was known that the signal component of the data 
was very large by comparison with the noise component. 

7 

The signal component is of the order of magnitude of 10 
meters; the noise of the order of magnitude of 1 meter. 

It was also known that the signal component is concentrated 
at the low end of the frequency spectrum. Based on these 
two facts, a linear smoothing, or spectral window, as 
described in Appendix A, Section A. 2.1, must be applied 
to the spectrum to prevent spectral leakage of the signal 
component into the higher frequency spectral bands. 

A series of test runs on actual range data, similar 
to those performed on simulated data and discussed in 
Section A. 2. 2, were performed to determine the amount of 
smoothing required. The linear smoothing selected was the 
application of the hanning window, three times. For 
these runs quadratic smoothing, by averaging of spectral 
estimates over a frequency band, was not used because the 
highest possible resolution on the spectrum was desired, 
at the expense of a large variance in the spectral esti- 
mates of the noise content. (Appendix A, Section A. 2.1 
and Lecture 5 of Reference 1 present more detail on the 
effects of linear and quadratic smoothing of spectral 
estimates) . 
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Data runs were then made using the selected smoothing 
method on data from 23 satellite passes, to determine the 
effective signal cut-off frequency. Figures 1, 2, and 3 
present typical results of the computations. From each 
satellite pass, 3000 data points at 10/second were used 
in the spectral estimate, yielding a frequency resolution 
of .00333 hertz. Only the first .333 hertz of the total 
frequency range (0-5 hertz) is plotted, as the region 
of interest in the determination of signal cut-off is at 
low frequencies. Figure 1 presents range data, Figure 
2 elevation data and Figure 3 azimuth data from GEOS 
Rev 5050, taken by the Wallops AN/FPQ-6 radar (NWALI3). 

The spectral density estimate for the range data, 
plotted in Figure 1 shows a well defined signal peak at 
0 frequency, about 125 db. The signal peak drops off 
smoothly, and at about .05 Hz merges with an apparent white 
noise at a level of -20 to -30 db . The spectral estimates 
of the noise are what would be expected from the the method 
used, which is not designed for optimum estimation of noise 
characteristics. The noise data scatter, almost 30 db, 
is in agreement with theory. The apparent short period 
trends are evidence of the serial correlation caused by 
application of the linear window function. 

All of the spectral density estimates for range data, 
performed using the above methods, showed the same charact- 
eristics. The cut-off frequency varied between .033 hertz 
and .05 hertz, and the average level of the noise spectrum 
varied at most 3 db, among the data taken from beacon track 
records. Spectral estimates performed on range data from 
a skin track mission, using the same computational methods 
show similar results, as shown in Figure 4. The noise level 
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AZIMUTH DATA, GEOS REV 5050, NWALI3 
SIGNAL CUT-OFF DETERMINATION 
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is about 30 db higher than that of the beacon track, and 
the signal cut-off frequency is about .033 hertz. 

Spectral density estimates of the elevation data 
and azimuth data, as shown in Figures 2 and 3 show a sig- 
nal cut-off at about .033 hertz and a noise level of 
-90 to -110 db. These results are also typical of 
all of the data examined. 

The conclusion to be drawn from this study is that 
the signal content of the data is effectively contained 
in frequencies below .05 hertz. It is possible, and 
indeed probable, that this conclusion is conservative. 
Since the windowing or linear smoothing method used 
has the effect of broadening any spectral peak, the 
signal cut-off frequency can be as much as .01 hertz less 
than that shown. 


2.1.2 Beacon Track Spectral Analysis 

The studies described above led to adoption of a 
standard method for obtaining power spectral density es- 
timates from the radar measurement data. The method is 
as follows: 

a) 3000 data points, sampled at 10/second, were 
used for all spectral estimates. The nyquist, 
or folding frequency, is 5 hertz for data sam- 
pled at 10/second, and the frequency resolution 
of the raw periodogram is ,00333 hertz. 

b) The data was not prefiltered to remove the signal 
content. Three successive applications of the 
hanning frequency window were done, to decrease 
spectral leakage from the signal component of 
the data. 
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c) The raw periodogram was smoothed to reduce noise 
variance by averaging over 15 adjacent points. 

This reduced the spectral resolution (interval 
between data points on the frequency axis) to 
.05 hertz. The number of degrees of freedom 
of the estimate of each spectral point is 30, 
and the spread, in db, including 801 of the 
data points is about 3 db. (See Reference 1, 
page 5-18) . 

Data from 20 revolutions of GEOS-II, taken from 
measurements by the Wallops AN/FPQ-6 and AN/FPS-16 radars, 
were reduced by this method. The results for GEOS rev 
5050 are plotted on Figure 5 through 10. 

The power spectral density of the range measurements 
on rev 5050, from the Wallops AN/FPQ-6, are plotted on 
Figure 5. The plot shows the signal spectral peak at 0 
db and the noise spectrum essentially flat from .1 hz to 
5.0 hz. This spectrum is typical of all of those examined. 
The similar spectral density estimate for range data from 
the AN/FPS-16 is shown on Figure 8. 

Figure 6 shows the spectral density estimate for 
the elevation data from the Wallops AN/FPQ-6. This plots 
shows a noise spectrum falling off linearly about 10 db 
between .05 hertz and 2 hertz, then decreasing another 
3 db out to 5 hertz. Two frequency spikes appear at exactly 
2 and exactly 4 hertz. It has been determined that these 
spikes are caused by an eccentricity in the elevation 
encoder mechanism, which creates errors with a fundamental 
frequency of 8 hertz, plus harmonics at 16, 24 and 32 hertz. 
The 8 and 32 hertz components fold to 2 hertz, and the 
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16 and 24 hertz components fold to 4 hertz in the sampled 
data. The spectral density estimate for elevation data from 
the Wallops AN/FPS-16, shown on Figure 9, is similar in 
general nature, but does not show the frequency spikes. 

The spectrum of the azimuth data from the Wallops AN/FPQ-6 
is shown on Figure 7, and that from the AN/FPS-16 on Figure 
10 . 


2.1.3 Skin Track Data Spectral Analysis 

Since essentially all of the GEOS-II data was taken 
in beacon track, no skin track data from the GEOS-II was 
analyzed. Data was selected from skin track measurements 
by the Wallops AN/FPQ-6 on Object 1314, and the spectral 
estimates were obtained by the method described in the 
previous section. 

A plot of the spectral density of the skin track range 
data is shown on Figure II. This spectrum shows a pronoun - 
ed peak in the noise spectrum at .2 hertz. The noise level 
drops off about 20 db between .2 hertz and 1.2 hertz, then 
remains essentially level out to 5 hertz. It is assumed 
that the object, which is a large booster motor, might be 
tumbling or precessing with a period of about 5 seconds, 
causing a peak at about .2 hertz in the noise spectrum. 


2.2 DATA COMPACTION STUDIES 

The data compaction studies were carried out using 
the COMPACT program which is described in Appendix A. The 
goal of these studies was to investigate the properties 
of the filtered and compacted data, to better determine 
an optimum method for data compaction. 
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2.2.1 GEOS - II Rev No. 5462 


In this test computation, data from the Wallops 
AN/FPQ-6 on Rev 5462 was filtered but not compacted. The 
data was kept at 10/sec. for analysis purposes. The filter 
used was a 401 point (40 second) non-recursive low-pass 
filter with nominal cut-off at .15 Hz. The number of data 
points, N, the cut-off frequency, f c , and the sampling in- 
terval T are related by the equation 



given in Section A. 3.1. The cut-off frequency f = .15 
hertz is not optimum from the sense of removing all noise 
from the data. The spectral analysis studies reported 
above indicate that f c = .05, with N = 1201, would be the 
best choice, if it is desired to eliminate essentially 
all noise while retaining the signal intact. However, 

N - 1201 means the loss of two minutes of data, one min- 
ute at the start and one minute at the end of the data 
span. It was recommended by Wallops personnel that a 
higher cut-off frequency be selected to reduce the amount 
of data lost. The question of optimum interval selection 
is discussed further in Section 3.2. 

In order to obtain a comparison between filtered and 
unfiltered data, the unfiltered range data was used in an 
A/OMEGA run to determine an orbit and obtain residuals. 

For this run, the range RMS was .75 meters. The smoothed 
range data were then used through one iteration of A/OMEGA, 
with the orbit elements determined from the unsmoothed data 
to obtain residuals on the smoothed data. The RMS of resid 
uals of the smoothed range data from the orbit obtained 
from the unsmoothed data was .43 meters. 
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A complete plot of the approximately 9000 residuals 
was not made. 

Figure 12 is a plot of the residuals of the smoothed 
(filtered) range data with respect to the reference orbit, 
on a compressed time scale. Every 50th data point (one 
every 5 seconds) was plotted on Figure 12. Figure 12 gives 
an idea of the noise left in the data after filtering and 
shows a low frequency trending of the data. Figure 13 
shows the corresponding unsmoothed residuals. 

Figure 14, 15 and 16 are plots of short segments 
of the filtered and unfiltered residuals vs. time, plotted 
10/second on an expanded time scale. Figure 14 shows data 
just after acquisition, Figure 15, near PCA, and Figure 16 
near end of track. 

On these plots, the filtered data shows jumps of 
almost half a meter at the two ends of the tracking span, 
but is very smooth near PCA. Examination of the A/OMEGA 
run shows that the jumps appear when the seconds column 
of the time tag changes by one in the fourth decimal place. 
Also, the size of the jumps is almost exactly one one- 
thousandth of the range difference taken over one-tenth of 
a second. The changes in the time tags are caused by 
transit time corrections, and the time is carried on the 
input tape to only four decimal places. The jump occurs 
because the computed value of the range is computed at 
a time which is in error in the fifth decimal place, com- 
pared to the time at which the measurement was made. 
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These j:umps result in an error which is locally non- 
random and has a variance of one- twelfth of the mean 
square value of the jumps. This error is probably neg- 
ligible when compared to the other errors in the data. 

See, for example. Figure IS, and note that the noise 
on the data seems to be about as big near PCA, where no 
jumps occur, as it does near the ends of the plot, where 
the jumps are the largest. 

The jumps can be prevented in one of two ways. The 
first is to change the data format to carry five decimals 
on the seconds of time. The second is to compute the 
transit time correction in A/OMEGA, in double precision, 
and input non- transit-time-corrected time. 

The plot of residuals of filtered data with respect 
to the reference orbit, in Figure 12, shows a systematic 
or low frequency component with a period of about 7 minutes 
and amplitude of .5 meter. In order to determine whether 
this might represent some actual satellite motion, rather 
than a low frequency radar noise, data was reduced from the 
Wallops AN/FPS-16 on the same pass, using the same method. 
The residuals of the filtered range data using a reference 
orbit obtained from the unfiltered range data are shown on 
Figure 17. The low frequency component which appears in 
the AN/FPQ-6 data is not apparent in the AN/FPS-16 data. 

The tentative conclusion is that it is caused by low fre- 
quency noise in the radar range system. 

In the course of examining the residuals of the 
unfiltered range data from the orbit; plotted in Figures 
14, 15, and 16, an apparent short term systematic pattern 
was observed. 
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FIGURE 17 

RESIDUALS OF FILTERED RANGE DATA 
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The same pattern appeared in the residuals obtained 
by differencing the unfiltered data and the filtered range 
data. Figures 18, 19 and 20 show the distinctive pattern 
of the differences between filtered and unfiltered range 
data, when plotted vs time. Instead of plotting 9000 
data points, three ten second intervals (100 points each) 
were plotted as typical examples. One set near acquisition, 
one near PCA and one near loss of track, was plotted. All 
three plots show points with a nearly linear distribution be- 
tween plus and minus one meter, with sporadic wild points 
outside those limits. In addition, the first and third of 
the three show a pattern in the successive data points. 

The linear distribution and the patterns shown can 
both be explained by assuming that the predominant source 
of error in the range data is granularity caused by the 
low order bit size. This granularity or truncation error 
would be linearly distributed about the mean. 

When the range rate is almost constant, as it is 
at the two ends of the track, the successive errors will 
follow a pattern of a saw-tooth curve with slope depending 
on the velocity and slowly changing as the R- changes. 

When the R is exactly a multiple of the granularity, 
the errors will become nearly constant for a few points; 
then as the R changes, the sawtooth pattern resumes. 

Near PCA, the range acceleration is great enough that 
these patterns do not form. 

As a test of the explanation, the effect of data granu- 
larity was simulated by computing the exact range to the com- 
puted orbit, at 10/second over a 6 second period, then intro- 
ducing the granularity caused by a low order bit weight of 
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Time (Minutes). 




1.786 meters. The residuals obtained are plotted on Figure 
21, along with the actual residuals from the orbit. 

The similarity in the short period trending in the 
two curves is apparent. The differences, for example 
in the interval from 22.4 to 23.0 seconds, are caused by 
a bias i n the level at which granularity occurs. If one 
bit weight is added to the simulated residuals in that 
interval, the curves will be similar. 


2*2.2 Data Analysis on WICE Test 80 

In order to further investigate the possibility 
of detection of orbital variations with periods of the 
order of a few minutes, another reduction was made, us- 
ing data from both the AN/FPQ-6 and the AN/FPS-16 radars 
at Wallops Island, Range data from both radars, tracking 
WICE Test number 80, were filtered and compacted using 
the COMPACT program. The filter used was a 401 point 
(40 second) low pass filter with cut-off at 0.15 Hz. The 
filtered data were selected at every 5 seconds to elimin- 
ate serial correlation. 

The filtered data were used in an A/OMEGA short arc 
reduction, with the AN/FPQ-6 data used to determine the or- 
bit. Range residuals are plotted on Figure 22. The large 
range bias was also observed on the reduction of the normal 
WICE. test data. From inspection of the plotted data, no 
correlation can be observed between results from the two 
radars. The conclusion is that there are no systematic 
effects which can be attributed to short period orbital 
variations. 
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As a further comparison of filtered and unfiltered 
data, orbits were determined, using data from the AN/FPQ-6 
radar only. 

One orbit used the unfiltered range data, selected 
every five seconds. The other orbit used filtered range 
data, as described above, also selected every five sec- 
onds. The results are summarized in Table 1, which tabu- 
lates orbital elements from unfiltered data, filtered 
data, the differences of corresponding elements, and the 
estimated standard deviation of each element. The RMS 
of the residuals of range, azimuth and elevation for each 
orbit are also tabulated, along with the input a priori 
standard deviation of each variable. It is evident that 
: the two sets of elements are in agreement to well within 
the estimated error limits. 

2.2.3 NBER05 Data Comparisons 

Range, azimuth and elevation data from NBER05 on 
two passes, Rev 4910 and Rev 5000 were filtered and com- 
pacted using the COMPACT program. The filter used was a 
401 point (40 sec.) low pass filter with cut-off at .15 
Hz. The filtered data were selected at 1/5 second, and in- 
cluded in a three week long arc reduction. Normally 
processed data from the same passes were included in another 
reduction. Residuals from the two reductions are plotted 
on Figures 23 to 28. The residuals of the filtered data 
from the reference orbit show essentially the same trends, 
but significantly less short period variation. 
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COMPARISON OF ORBITAL ELEMENTS DERIVED FROM FILTERED 
AND FROM UNFILTERED DATA FROM NWALI3 ON WICE TEST 80 



UNFILTERED 

DATA 

FILTERED 

DATA 

A 

a EST 

X CM) 

- I62408S . 0 

-1624101.5 

-16.5 

71. 

Y (M) 

- 7366614 i 4 

-7366606.8 

7.6 

30. 

Z (M) 

1691209.5 

1691206.0 

-3.5 

23. 

X(M/S) 

-1671.3504 

-1671.3260 

.0244 

.1368 

Y (M/S) 

1690.0173 

1690.0135 

.0038 

.0219 

Z (M/S) 

6766.9897 

6766.9991 

-.0091 

.0390 

RMS (R) 
meters 

.762 

.388 

° INPUT 

2;m 


RMS (A) 
sec of 
arc 

16.5 

16.4 

50.0 


RMS (E) 
sec of 
arc 

21.6 

20.0 

50.0 



TABLE 1. 
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2.2.4 Long Arc Comparison 

In order to examine the effects of use of filter- 
ing and compacting radar data on the orbital solution 
and the residuals obtained from the data, a long arc 
solution comparison was made. Data from seven GEOS-II 
passes over a period of slightly greater than 24 hours, 
as listed in Table 2, was reduced twice; 

once using filtered range data and once unfiltered range 
data. The filter used was a 601 point low pass filter, 
with cut-off frequency of .1 hertz. The filtered data 
was selected at 1 point per 10 seconds, then reduced in 
the nominal manner. The unfiltered data was also selected 
at 1 point per 10 seconds and reduced in the same manner. 

In order to evaluate the accuracy of the orbital 
elements obtained by the two methods, the elements were 
used to predict measured values, and compared against 
actual measurements. For each set of elements, a 10 day 
orbit was computed, and data from 17 stations, a total of 
9175 measurements, were used to compute residuals. 

The RMS of residuals from the orbit determined by 
the unfiltered data was 180.8 meters, and from the orbit 
determined by the filtered data was 159.4 meters. 

I 

The indication is that the two orbit determinations 
are essentially the same. A true statistical comparison 
is difficult, because the effects of unmodelled errors 
prevents a valid estimate of what the standard deviation 
should be. 

The orbital elements and the RMS residuals values 
are listed in Table 3. 
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Residuals of the unfiltered data are plotted on 
Figure 29, and those of the filtered data on Figure 30. 


DATA USED IN FILTERED -UNFILTERED LONG ARC 


STATION 

REV 

DATE 

START TIME 

STOP TIME 

NBER34 

4960 

690201 

085950 

090620 

NWALI3 

4960 

690201 

090000 

090620 

NWALI3 

4961 

690201 

104700 

105420 

NBER34 

4961 

690201 

104650 

105420 

NWALI3 

4967 

690201 

224140 

224830 

NWAL34 

4967 

690201 

224150 

224920 

NWALI3 

4973 

690202 

092110 

092700 


TABLE 2 


COMPARISON OF ORBITAL ELEMENTS FROM FILTERED 
AND FROM UNFILTERED DATA ON LONG -ARC SOLUTION 



UNFILTERED 

DATA 

FILTERED 

DATA 

A 

X(M) 

-5541910.7 

-5541934.7 

-24.0 

Y 00 

-3939388.2 

-3939395.1 

6.9 

Z(M) 

3942974.9 

3942961.1 

13.8 

X(M/S) 

1331.5775 

1331.5589 

.0186 

Y (M/S) 

3668.2914 

3668.2816 

-.0003 

Z CM/S) 

5872.6921 

5872.6874 

.0047 

RMS {M) 

10.2 

14.3 


RMS 

PREDICTED 

180.8 

159.4 



TABLE 3 
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2.2.5 CSP Studies on Wice Test 75 

In connection with the CSP studies on WICE Test 
73, the range rate values from the CSP unit were filter- 
ed, but not compacted, before computing the low-noise 
CSP ranges. The filter used was a 1201 point low-pass 
filter with cut-off at .05 hertz. The result was a 
reduction in RMS of fit to a 20 second orbit, from 
7 . 6mm to 0 . 6mm. 

2.3 AUTO -COVARIANCE FUNCTION STUDIES 

Severval studies of the nature of the auto-covariance 
function of the residuals of unfiltered vs filtered range 
data were completed. The program used for these studies 
is subroutine FP9, (See Appendix A) which can be called as 
an option in the COMPACT program. FP9 produces a periodo- 
gram of a set of data, then Fourier transforms the periodo- 
gram to obtain the auto-covariance function of the data. 


2.3.1 Auto- covariance Function of Normal Beacon Data 

Data from NBER34 on Rev 4960 was processed by the 
COMPACT program, using a 601 point filter with cut-off 
at .1 hertz. The auto- covariance option was used. The 
periodogram computed was essentially featureless. A 
plot is not shown, because the periodogram contains over 
1000 points and is a random scatter of points, above .1 
hertz. Below .1 hertz the periodogram is essentially 
zero, since the residuals from wich the periodogram was 
taken have the frequencies from 0 .1 hertz filtered out. 
The auto-covariance function shows a spike of .477 meters 
at zero lag, agreeing with the computed variance of the 
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residuals. All other points on the auto-covariance function 
are less than .04 in absolute value, and appear to be com- 
pletely random. 

2.3.2 Auto-covariance Results on The Ship Test Data 

Data from the VANGUARD radar on Puerto Rican Trench 
Ship Test number 18 was treated in the fashion described 
above. 2000 points of range data at 10 per second were 
passed through COMPACT program with the option of computing 
the periodogram and auto covariance function. It was sus- 
pected that the autophasing mechanism was introducing 
spikes in the range system at one second intervals. The 
residuals did not conclusively show the spikes, as the 
general noise level tended to obscure them. The periodo- 
gram of the data did show frequency spikes, occurring at 
exactly 1., 2., 3. and 4. hertz, of almost equal ampli- 
tudes. The variance attributable the four spikes is 
.26 mtr. 2 of a total variance of 4.78 mtr. It would appear 
that the error introduced by this effect is negligible 
with respect to the other contributing factors to the 
noise. The 1 per second spikes are quite apparent in the 
plot of the auto- covariance function of the residuals. 


2.3.3 Auto-covariance Results on PMRMR1 Data 

Power spectral density estimates on range data from 
PMRMR1 had indicated that some peculiarity would be found 
in the data. In addition, residuals from short arc data 
indicated periodic discontinuities in the data. Range 
data from Rev number 4989 were processed by the COMPACT 
program with option to compute the auto-covariance and 
periodogram. 
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The resulting residual plots showed discontinuities 
occurring approximately every 11.9 seconds. Figure 31 shows 
a 20 second sample of the residual plot. The periodogram, 
the low frequency part of which is plotted on Figure 32, 
shows spikes at the harmonics of this fundamental period. 

The fundamental, at .084 hertz, is relatively weak because 
of the action of the filter, which acts like a .1 hertz 
high-pass filter on the residuals. The auto- covariance 
function, plotted on Figure 33, clearly shows the periods 
of the fundamental and first harmonic. 


2.3.4 Auto-covariance Results on Residuals from Short-Arc 

Data 

A program, AUTOCORR, (See Appendix A) was written 
to read data from the A/ OMEGA Plot output tape, and to apply 
FP9 to obtain the auto-covariance function and periodogram 
of the resulting data. 

Raw range data at 10/second from the Wallops AN/FPQ-6 
radar was reduced in a short-arc reduction and the resulting 
output tape representing residuals was used as input to 
AUTOCORR. Figures 34, 35 and 36 show the results. 

The raw range data was then filtered by a 601 point 
low-pass filter, with cut-off at .1 hertz, and another 
short-arc reduction performed. The results of using the 
residuals of the filtered data as input to AUTOCORR are 
shown in Figures 37, 38 and 39. 

The residuals of the unfiltered data, plotted on 34, 
show noise plus low frequency trending. The residuals of 
the filtered data, plotted on 37, show only the trending. 
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RANGE RESIDUAL VS TIME 
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The power spectrum (periodogram) of the residuals 
from unfiltered range data, Figure 35, shows the low fre- 
quency content cutting off at about .2 hertz, with white 
noise at higher frequencies. The periodogram of the resid- 
uals of filtered data. Figure 38, show about the same low 
frequency content but also show the white noise cut-off 
at about .1 hertz, caused by the filtering. 

The auto- covariance functions of the two sets of 
data show a very similar overall structure. Both show a 
first zero crossing of about 25 seconds, which is compati- 
ble with a signal cut-off at about .02 hertz. 


3.0 RESULTS AND CONCLUSIONS 


3.1 RESULTS 

3.1.1 Spectral Analysis 

The major results of the investigation of the 
spectral content of the radar measurement channels are 
as follows: 

a) The spectral density of the signal component 
of the data from range, azimuth and elevation 
measurement channels, tracking a satellite, 
is concentrated at frequencies lower than .05 
hertz. A low-pass filter, with sharp cut-off 
at .05 hertz, will retain all signal content 
of the data, and reject almost all noise con- 
tent. 

b) The spectral density of the noise content of 
the range measurement channel is essentially 
flat from .05 hertz to 5 hertz, for 10/second 
data. Most of this noise content is caused 
by signal granularity, or low order bit trun- 
cation. 

c) The spectral density of the noise content of the 
azimuth and elevation measurement channels are 
similar showing relatively flat noise to 1 hertz, 
and a decrease of about 10 db to 3 hertz, then 
flat noise to 5 hertz. 
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d) The spectral density from all satellite passes 
examined were consistent with each other. Var- 
iation from the described pattern always indi- 
cated bad data, either due to processing errors 
or radar malfunction. 


3.1.2 Data Compaction 

The major results of the investigation of techniques 
for data compaction of the radar data on GEOS tracking 
records are as follows: 

a) Filtering before selection of data improves the 
RMS of the residuals to a short-arc fit, over 
simply selection of raw data. The example given 
in Table 1 is typical. In that test, data from 
WICE test 80, the RMS of the residuals from the 
unfiltered range data was .762 meters, and from 
the filtered range data, .388 meters. This im- 
plies that noise with an RMS value of .656 meters 
has been filtered out of the data. An RMS 
of .622 meters for the differences be- 

tween filtered and unfiltered data was computed 
by the COMPACT program when the data was filter- 
ed. 


b) The improvement of RMS is a non-linear function of 
the filter cut-off frequency. Because, as was 
shown by the spectral analysis studies, the power 
spectral density of the noise in the range channel 
is flat above .05 hertz, the pure noise variance 
will be directly proportional to the filter cut- 
off frequency. The variance of the unfiltered 
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data in the example is .2323^ = ,0540m? The 

2 

variance of the filtered data is .1183m = .0140 

m . The remaining variance, .0400 m is 
attributable noise uniformly distributed between 
.05 and 5 hertz. Figure 40 is a plot of the 
RMS of fit, vs filter cut-off frequency, based 
on these assumptions. 

c) Filtering the data by non-recursive low-pass 
filters causes a data loss, on each end of the 
data span, inversely proportional to the filter 
cut-off frequency. Figure 41 is a plot of data 
loss vs cut-off frequency for 10/second data. 

d) Data filtered before selection produces essen- 
tially the same orbital elements as data select- 
ed without filtering, for both short-arc and long- 
arc orbital reduction. This demonstrates that 
filtering does not change the statistical content 
of the data, and that the accuracy of the data is 
unaffected by filtering. 

e) Short-arc orbital reductions using unfiltered range 
data show distinctive short period trends in range 
residuals caused by interaction of range rate with 
the data granularity caused by low order bit size. 

f) Short- arc orbital reductions using filtered range 
data show low frequency trends, in the frequency 
range 0 to .025 hertz. Comparison of data on 

the same satellite pass, taken by adjacent radars, 
indicates that the low frequency trending in the 
residuals from different radars are dissimilar. 
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This implies that the trends are caused by low 
..frequency radar noise, not by orbit modelling 
errors or by data processing methods. 

3.1.3 Serial Correlation and Frequency Folding 

a) On the basis of theory, white noise with a cutoff 
frequency of f is uncorrelated when sampled at 
1/ ( 2 f . On this basis, if the noise spectra of 
the radar data is essentially flat to 5 hertz, the 
noise is uncorrelated when sampled at 10/second. 

However, if the low frequency trending shown in the 
short-arc residuals is radar noise, then the total noise 
spectra of the radar cannot be considered to be 

flat, but is peaked at low frequencies. 

b) Computation of the auto covariance function of 
residuals from short-arcs, indicates that zero 
correlation occurs at a lag of about 20 seconds, 
which is probably due to low frequency noise 
with cut-off at about .02 hertz. 

c) Serial correlation of residuals is affected 
very little by fitlering of data. On the test 
described in paragraph 3.3.4, the lag for neg- 
ligible serial correlation is approximately 
twenty seconds for both filtered and unfiltered 
data . 

d) On unfiltered data, any sampling rate less than 
10/second causes frequency folding. On filtered 
data, the sampling rate must be less than 1/(2 f c ) 
to prevent frequency folding. 

e) The theoretical optimum sampling rate to avoid both 
frequency folding and serial correlation is 1/(2 f C D * 

This is applicable only to data with white noise 
having a cut-off frequency f c . 
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3.2 RE COMMENDAT I ONS 


3.2.1 Data Compaction Method 

a) Based on the above results, it is recommended that 
the radar data which is to be stored for general 
scientific use be filtered with a low pass filter, 
with cut-off frequency at .1 hertz, then sampled 
at one sample per ten seconds. The .1 hertz fil- 
ter is recommended to reduce the random noise 
content essentially to zero, and to provide for 
negligible error contribution caused by frequency 
folding or aliasing when the data is sampled at 
one per ten seconds. With this filter, 30 seconds 
of data will be lost at each end of each segment of 
data. This is considered acceptable in terms of 
total data content. Any lower frequency cut-off, 
with corresponding greater data loss, would not be 
acceptable. The one per ten second sampling is 
recommended as the largest interval which simultan- 
eously 

1) Provides negligible serial correlation for 
short-arc residuals 

2) Provides negligible effects of frequency fold- 
ing for a .1 hertz cut-off filter. 

b) Data which is to be used only for long-arc reduct- 
ions can be compacted by simply selecting data at 
1 point per ten seconds or 1 point per 20 seconds. 
The slight increase in bias and low frequency 
error which can theoretically be expected to result 
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from frequency folding, apparently does not 
appreciably affect the orbital reduction. 

c) Data should not be sampled at greater than 

20 second intervals, pending further investi- 
gation, as the frequency folding of the low 
frequency noise might appreciably affect the 
data . 


3.2.2 Spectral Analysis 

a) The consistency of the spectral estimates is so 
good that it is recommended that a regular, 
standard check on radar operation should be an 
examination of the spectral density curves of the 
data. This check could be performed on raw 
data, before any reduction is done, and if done 
as a standard operational procedure would serve 
to identify instrument malfunctions early in the 
operation. 
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APPENDIX A 


COMPUTER PROGRAMS DEVELOPED 



A. 0 COMPUTER PROGRAMS DEVELOPED 
A. 1 SUBROUTINE PACKAGE 

In order to implement the methods developed for Spec- 
tral Analysis and data compaction, a number of computer 
programs and associated subroutines were programmed. Since 
many of the numerical methods used are applicable to more 
than one portion of the study, and in fact are generally 
applicable to many other data processing problems, a set of 
general purpose subroutines were programmed to implement 
these methods. The subroutines are tailored to emphasize 
the methods which are most useful in analysis of the C-Band 
radar data, but are written so as to be generally applicable 
to many frequency analysis problems. Reference 2 describes 
the subroutine package in detail and Reference 1 discusses 
the theory involved. 

A. 1.1 Subroutine Functional Descriptions 

Five of the subroutines are used to apply the, fast 
Fourier transform algorithm to data arrays. The finite 
Fourier transofrm, implemented by the fast Fourier trans- 
form algorithm, is a component of most programs for power 
spectral density estimation, digital filtering and auto- 
covariance estimation.. The subroutines are as follows: 

FP 6 Factors a composite (non-prime) number into 
prime factors and computes partial products 
of the factors for use in the fast Fourier 
transform (fFT)..- The fFt algorithm used is 
not limited to a number of points which is 
a power of two.. It is applicable to any 
composite number for which the prime factors 
are available. 
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FP 7 Computes the fast Fourier transform of a 
complex data array. Either the direct or 
the inverse transform can be computed. 

FP 8 Re-orders the output of FP 7, stored in 
digit inverted order, into normal order. 

FP10 Computes the array of complex exponentials 
used in FP 7. 

FP15 Used to economize on computing time and array 
storage, by allowing the transform of real 
arrays "packed” into complex arrays. FP15 un- 
scrambles the transform of the packed array to 
produce the transform of the real array. 

These five subroutines, called in the appropriate sequence, 
will perform any desired finite Fourier transform operation. 

FP 1 Computes the coefficients of a non^recursive 
digital filter over N points, to obtain a 
high-pass, low-pass, band-pass or band stop 
filter . 

FP 2 Applies the Lanczos sigma smoothing as a lag 
window function to a non-recursive digital 
filter array, to minimize the Gibb’s effect 
oscillations . 

FP 3 Computes the frequency response function 
(transfer function) of a symmetric non- 
recursive digital filter. 



FP 4 Applies a non-recursive digital filter to a 
data array. 

FP. 5 Modifies the coefficients of a non-recursive 
digital filter so that the gain of the filter 
is unity at a specified frequency and zero 
at another specified frequency. 

FP18 Applies a non- recursive digital filter to 

a data array by using the finite Fourier con- 
volution property. The filter and the data 
to be filtered are Fourier transformed and 
the product of the transforms is inverse 
transformed to obtain a filtered array. 

Subroutine FP11 applies a recursive three point notch filter 
to a data array, to remove isolated frequency spikes from 
the; spectrum. I 

Four subroutines are used to generate test data, of known 
frequency content and statistical properties, to provide 
test cases for the other subroutines. 

FP12 Computes a sequence of pseudo-random numbers 
from a uniform distribution, zero to one, 
by the multiplicative congruence method. 

FP13 Filters a sequence of pseudo- random numbers 
to obtain pseudo-white noise with specified 
amplitude and cutoff frequency. 

FP14 Constructs a data array to simulate s signal 

with up to ten frequency spikes (pure sinusoids) , 
of specified amplitudes and frequencies, added 
to white noise. 
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TDATA Computes simulated radar tracking data (range, 
azimuth, and elevation) from a simplified or- 
bital tracking model to provide idealized 
error free, data inputs to the subroutines. 

Since the theory underlying all of the methods used is based 
on equally spaced data, and since the data contains wild 
points which must be edited out, an editing subroutine is 
used. The editing subroutine is FP16, which takes the 
Fourier transform of the array, deletes high frequency com- 
ponents, then inverse transforms to obtain a smooth function. 
The smooth function is compared point by point with the 
original data, and all data points differing from the corre- 
sponding smoothed value by more than an error limit are 
edited and replaced by the smoothed value. 

Subroutine FP20 is used to subtract a double precision 
quadratic function from an array X, store the coefficients 
of the quadratic, then add the same quadratic to the 
modified array. This is used to remove the large magnitude 
data component from data before filtering, and to restore it 
after filtering, in order to minimize the effects of trun- 
cation error in the computation. 

Subroutine FP 9 computes the periodogram and autocovariance 
function of a real data array. 

Subroutine FP19 computes the cross spectrum and cross- 
covariance of two real arrays. 
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Subroutine FP17 computes an estimate of the power spectral 
density of a real data array. The spectral estimate is 
obtained by computing the discrete Fourier transform of 

the dat^a, then averaging and squaring to obtain smoothed 

\ 

spectral estimates. 

Three\ types of averaging are used: 

1) The data is segmented, before transforming, into 

KSEG aubarrays of NSEG points each. After transforming 
and squaring to obtain raw periodograms , the spectral 
estimates (periodogram values) from corresponding fre- 
quency points are averaged. 

2) The Fourier transforms of the data segments are averaged, 
before squaring and averaging, by applying the hanning 
spectral window KWIND times. 

3) The squared estimates from each segment are averaged over 
KAVG points, and only the midpoint of each KAVG points 

is stored. 

A. 2 SPECDEN PROGRAM 

The SPECDEN program was developed to compute estimates 
of power spectral density of radar measurement data. 

A. 2.1 Fucntional Description 

The SPECDEN program either accepts as inputs radar 
measurement data in the Wallops PASS1 tape format or on an 
input option, generates test data internally for the purpose 
of program testing. If the PASS1 tape is used as input, 
any word from the data record can be read into storage for 
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of estimating the power spectral density. The start and 
stop times of the data input are also specified on input • 
option. If sporadic gaps occur in the input data, zeroes 
are substituted into the stored array. An editing sub- 
routine (FP16) is used to substitute smoothed data points 
in for these missing points, and to edit out wild points 
and substitute smoothed data in their place. 

Three types of internally generated test data are avail- 
able, on input option. 

a) Simulated noise free radar measurement data, from a 
simplified satellite tracking model. 

b) Up to ten sine and cosine functions with amplitudes 
and frequencies as input options. 

c) White noise with a specified input r.m.s. value, 
filtered to a specified cutoff frequency. 

These three input options can be used in any combination. 

The data, once read in (or generated internally) and edit- 
ed, can be filtered, if desired, on input option, as follows: 

a) Up to 12 recursive notch filters with specified 
center frequencies and bandwidth can be applied 
to the data, using FP11. 

b) A non-recursive high-pass, low-pass or band-pass 
filter can be applied to the data using FP18. 
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These filters can be used to remove known frequency compon- 
ents from the data before estimating the spectral density. 

The spectral density estimation is done by FP17. This 
subrouitne obtains the spectral density by taking the finite 
Fourier transform of the data to obtain the spectrum of the 
data, then squaring each spectral component to obtain the 
periodogram of the data. Three types of spectral averaging, 
or smoothing can be applied to the data to obtain an esti- 
mate of power spectral density. 

a) Linear smoothing of the Fourier transform of the 
data, performed before squaring or averaging, ob- 
tained by applying the hanning spectral window. 

This smoothing improves the resolution of the 
spectral estimate by decreasing spectral leakage 
from non harmonic components. 

b) Quadratic smoothing of the periodogram, performed 
after squaring, by averaging a number of adjacent 
values of the periodogram. This quadratic smooth- 
ing improves the variance of the spectral estimates 
of random noise. 

c) Ensemble averaging of the periodogram, obtained by 
segmenting the. data, transforming each segment and 
obtaining the periodograms , the averaging corres- 
ponding frequency components. 

The effects of these three types of smoothing on sample 
data are discussed in the next section. 

The theory of linear and quadratic smoothing of spectral 
estimates is discussedin Lecture 5, Reference 1. 
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The output of SPECDEN is a listing and printer plot 
of the estimated power spectral density of the input data. 


A. 2. 2 Test Data Output 

The data presented in Figures A-l through A- 7 demon- 
strates the use of linear smoothing of the periodogram to 
obtain spectral estimates. For each test 200 data points 
were generated at a simulated time increment of 1/second. 

The result is 101 spectral points from 0 to 0.5 hz, at a 
spacing of .005 hz. All the periodogram from 0 to 0.5 hz 
is symmetric with that from 0 to 0.5 hz, and in accordance 
with convention, the total power is presented in the positive 
portion of the spectrum. 

All spectral data are presented in decibels, obtained 
by taking 10 log 10 of the spectral points. 

Test 1 data consisted of a steady state component 
and periodic components at .1, .2, .3, .4 and .5 hz, all 
with amplitude 1. Figure 1 shows the resulting spectra, 
with spikes at the input frequencies at 0. db for the 
steady state component and -3 db for the periodic components. 
This test demonstrated that the program performs correctly 
on idealized periodic, harmonic data. 

The computational noise level, caused by truncation 
in the 27 bit mantissa of the floating point word in the 
computer, is about -140 db, so the dynamic range of the 
spectral estimates if 140 db. 
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FIGURE A- 
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Test 2, with results shown in Figure 2, used the 

same periodic components, with a steady state component 
7 

of 10 . The truncation error increased to about -50 
db, the steady state spike increased to 140 db, as it 
should, and the periodic component spikes stayed constant 
at -3 db. 

Test 3 presents a graphic demonstration of spectral 
leakage caused by a large non-harmonic component of data. 
The input was a component with amplitude of .0025 hz, or 
one half of the fundamental frequency of the data interval, 
and components with amplitude 1 at .1, .2, .3, .4 and .5 
hz, as before. The spectral leakage caused by the large 
non-harmonic component completely hides all other detail in 
the spectral estimate. 

Figures 4, 5, 6 and 7 show the effect of increasing 
amounts of linear smoothing on the spectrum of the same 
data as presented in Figure 3. In Figure 4 the hanning 
window was applied to the transformed data once before 
squaring. In Figure 5 the hanning window was applied 
twice, in Figure 6, four times, and in Figure 7, eight 
times. The effects to be noticed are these: 

a) As the hanning window is applied more times, 
the peak near zero hertz gets sharper until 
Figure 7, where it is broader. It appears 
that about four applications of the hanning 
window is optimum, for improving resolution 
of peaks of this magnitude. 
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b) The signals at .1, .2, .3, .4 and .5 hertz emerge 
from the frequency leakage, but are broadened and 
decreased in amplitude by the effects of the window. 

The data presented in Figures A-8 through A-14 pre- 
sent the effects of quadratic smoothing on spectral density 
estimates computed from random noise. The input data for 
all of the these results were computed by generating a set 
of random numbers with a uniform distribution, then applying 
a linear filter with cutoff at .25 hz. In addition, a fre- 
quency spike, consisting of a cosine function with frequency 
0.1 hz and amplitude 1 was added to the data. A data samp- 
ling rate of 1/second was used, resulting in a frequency 
range from 0. to 0.5 hz in the output spectral estimates. 

The power spectral density estimates shown were 
arrived at by computing the periodogram over N data points, 
then averaging over M frequency intervals. The numbers N 
and M are given below: 


Figure N M 


A-8 

100 

1 

A-9 

200 • 

2 

A- 10 

500 

5 

A- 11 

1000 

10 

A- 12 

2000 

20 

A- 13 

5 000 

50 

A-14 

10000 

100 
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FREQUENCY (HZ) 








The result, in each case, is 51 spectral points, 
with a frequency interval of .01 hz. The decrease in 
spread of the points with increased smoothing is evident. 
The data spread, or variance, according to theory should 
be inversely proportional to the number of points, M, 
used in smoothing. The results are in agreement with this 
The frequency spike at . 1 hz is not affected by the smooth 
ing, also in agreement with theory. 

These results were presented to demonstrate the 
effects of the linear and quadratic filtering used in 
obtaining the spectral density estimates and to demon- 
strate that the SPECDEN program performs correctly on 
artificially generated input data, for which the results 
are known from the theory. Other tests, not presented 
here, were conducted to validate the editing methods, 
the filtering methods, and the accuracy of the spectral 
density estimates on other combinations of artificially 
generated test data. 
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A. 3 COMPACT PROGRAM 

The COMPACT program has the function of preparing com- 
pacted radar data from the PASS1 raw data tapes. 


A. 3.1 Functional Description 

The COMPACT program accepts as input a magnetic 
tape in the Wallops PASS1 format as given in Table A-l. 
The output ;of the COMPACT program is a tape in the same 
format, with filtered and compacted data. The design 
philosophy was that the COMPACT Program would fit into 
the satellite data reduction operation without requiring 
any change in any of the other computer programs used. 
The program, and its use, satisfy this philosophy. If 
data compaction is desired, the PASS1 raw data tape is 
used as input to the COMPACT Program, and an output 
PASS1 filtered data tape is prepared, which is used as 
input to the ELPRO program. If data compaction is not 
desired, the PASS1 raw data tape is used as input to 
the ELPRO program. In order to accomplish this, the 
following operations are performed: 

a) The pre-mission calibration is written on the 
the output tape unchanged. 

b) All tracking data before an input start time 
and after an input stop time are deleted. 

c) The data is first filtered with a digital 
low pass filter. 
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Tape 

Format 

(output) 


Mode of 

Writing: Binary Density: 800 BPI 


Words/Record: 34 


Explanation 

of Record: 


Each 

word consists of 36 bits: 



1 

I 

Run Number 



2 

A 

Radar Number (0-11); Model Number (12-35) 



3 

A 

Radar Mode (0-5) 



4 

A 

Control Mode 


5 $ 

6 

DP 

Time (hours) 


7 § 

8 

DP 

Azimuth (degrees) 


9>§ 

10 

DP 

Elevation (degrees) 


11 U 

12 

DP 

Slant Range (feet) 


13 § 

14 

DP 

Delta Azimuth (degrees) 


15 § 

16 

DP 

Delta Elevation (degrees) 


17 § 

18 

DP 

Delta Slant Range (feet) 



19 

I 

Azimuth Error (bits or Temp bits) 



20 

I 

Elevation Error (bits) or Humidity (bits) 



21 

I 

Range Rate (feet per sec x 10000) 



22 

I 

Roll (bits) or Velocity of Propagation (feet per 

sec) 


23 

I 

Pitch (bits) or Doppler Freq Counter (Megahertz) 



24 

I 

Yaw (bits) or Transmitter Freq (Hertz) 



25 

I 

Ship Velocity (knots) or Temp A 



26 

I ■ 

Ship Latitude (degrees x 1000) or Humidity A 



27 

I 

Ship Longitude (degrees x 1000) 



28 

1 

AGC (bits) or Calibration- Number (bits) 


29 § 

30 

DP 

Elapsed Time (hours) 



31 

0 

Milgo Aux Bits Recorded Only (0-32); R bits (33-35) 


32 

0 

FPQ-6 Status Bits (0-29); Milgo Aux Bits Recorded 

and 


'* 


Transmitted Bits (31-35) 



33 

A 

Azimuth Lag Correction Bit (bit 5) or Calibration 

Data 


34 

A 

Elevation Lag Correction Bit (bit 5) 



I 

= Integer 



A 

= Alphanumeric 



DP 

= Double Precision 



0 

= Octal 



TABLE A- 1 

PASS1 RECORD FORMAT 
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d) Every NSEL^^ 1 point is then selected from the 
filtered data and written on an output tape. 

e) Each output record consists of a PASS1 record. 
All unfiltered (auxiliary) data from the PASS1 
are recorded unchanged. The items to be fil- 
tered are specified by input option, and any 
word on the record can be specified. Range 
only, or range, azimuth and elevation will 
usually be selected. 

f) The post-mission calibration is written on the 
output tape unchanged. 

The result is a simulated PASS1 tape, complete 
with calibrations, with the selected tracking data filtered 
and compacted. 

Several secondary functions are performed on the 
data, as follows: 

a) If a data glitch has caused the "model number" 
on pre-mission calibration to erroneously show 
post mission status, the model number is correct 
ed. 

b) The data is edited for wild points by testing 
the residuals of filtered data from the raw 
data, rejecting all data for which the residual 
exceeds an input multiple of sigma, and sub- 
stituting the filtered data point. 
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c) The power spectral densities of the raw data 
and of the filtered data are computed and 
output on a printer plot, for use as a check 
on data quality. 

d) The data residuals, filtered minus raw data, 
are computed and output as a printer plot. 

e) A histogram of the data residuals is computed 
and output as a printer plot. 

f) On an input option, FP9 is used to compute 

the periodogram and the autocovariance function 
of the residuals, and to output printer plots 
of the results. 

The data selected is filtered by a non-recursive 
low pass digital filter. The number of filter weights 
(NFILT) is input, as is the cutoff frequency (TH2) in 
hertz. The filtered output Y, corresponding to unfiltered 
data input X, with filter coefficients a, is given by 



where N = (NFILT-l)/2. 

Using this type of filter, N points at the beginning 
of the data and N points at the end are lost in the filtering 
process. A study was made to determine the number of filter 
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coefficients necessary for the low-pass filter functions 
to be performed. The following general conclusions were 
reached. The details of the theory are given in Lecture 
2, Reference 1. 

a) Every low-pass filter has a finite transition 
region on the plot of ..gain vs frequency, be- 
tween the pass-band, in which the gain is 
essentially unity and the stop-band, in which 
the gain is essentially zero. 

b) The width of the transition region is inversely 
proportional to the number of points in the 
filter. 

c) The cutoff frequency (the frequency at which 
the gain is 0.5) is the center of the trans- 
ition region. 

d) The number of points in the filter, NFILT, 
which determines the amount of data lost at 
each end, is inversely proportional to the 
cutoff frequency of the filter and to the 
data interval T. 


data interval T. 


NFILT = 6 . 0/ (T*TH2) +1 

Table A- 2 summarizes the relations for T = .1 
second, or 10/second data. 


Cutoff 

Frequency 

(Hz) 

Number 
of points 
in the filter 

- 

Data lost at 
each end of 
the filter 
(seconds) 

1 . 

61 

3 

. 5 

121 

6 

.25 

241 

12 

.15 

401 

20 

.1 

601 

30 

.05 

1201 

60 

.0167 

3601 

180 


TABLE A - 2 


RELATIONSHIP BETWEEN CUT-OFF FREQUENCY, 
LENGTH OF FILTER AND DATA LOST 
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